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Abstract 

In this article, we follow the trajectory of a photon in the human 
cornea. Prior to experimental evidence of the cornea's gradient refractive 
index (GRIN) nature, schematic eye models used a constant or average 
as the corneal refractive index. A few recent models respect the intra- 
corneal GRIN, but are based on the paraxial approximation, and thereby 
have limited validity in peripheral visual field analysis and non-paraxial 
photon tracking. Here, we introduce the Trajectron algorithm. It uses a 
closed form solution of the ray equation of constant axial GRIN media, and 
evaluates over piece-wise constant GRIN. Using Trajectron, we present 
the first quantitative demonstration of non-paraxial intra-corneal light- 
bending phenomena. This demonstration has significant implications in 
the development of refractive surgery algorithms. 
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1 Introduction 

1.1 Aim 

To determine the trajectory of a photon being observed in the human cornea, 
from a point within the anterior central corneal epithelium to a point within 
the posterior stroma. 

1.2 Motivation 

Refractive error places an enormous disease burden on society, with net esti- 
mates of over two billion people affected worldwide [H ITO l l45 l l48 l l47 l [T2 l [38 l 
HOl m [39l [H]. Myopia, hyperopia, presbyopia, astigmatism, and cataracts are 
among the most common forms of refractive error. Refractive surgery is an effec- 
tive and growing method of correcting refractive error. There are a high volume 
of cases performed annually, including cataract extractions with intraocular lens 
placement, laser-assisted in-situ keratomileusis (LASIK) by microkeratome inci- 
sion or by femtosecond laser incision, photorefractive keratectomy (PRK), and 



laser-assisted sub-epithelial keratomileusis (LASEK). Refractive surgery tech- 
niques and technologies have seen significant advances over the years including 
current reports of greater than 90% patient satisfaction post-operatively |42l |4] . 
However, outcomes remain suboptimal and patients are often still left to manage 
imperfect vision following surgery. Most modern refractive surgery techniques 
involve direct modification of the corneal epithelium and stroma. It is there- 
fore important to have a precise and detailed understanding of the refractive 
properties along this path. It is additionally imperative to understand how 
those properties govern the behavior of light both as a source of visual stimulus, 
and in amplified form, as lasers for refractive and other ophthalmic interven- 
tion. Furthermore, a precise quantitative mapping of the refractive properties 
of the human cornea is requisite for good design of keratoprosthetic implants. 
In this paper, we embark on a mathematical exploration of the photon tra- 
jectory from the anterior central corneal epithelium to the posterior stroma. 
We operate within the realm of geometric optics, which assumes the limit of 
short wavelength, where light propagation can be approximated by rays and 
studied using geometry. We introduce the Trajectron algorithm for comput- 
ing the photon trajectory in a model respectful of the cornea's axial GRIN 
nature. Meridional intra-corneal GRIN is of significantly less magnitude than 
axial intra-corneal GRIN, and is assumed to be zero in this model. Trajectron 
uses a closed form solution of the ray equation of constant axial GRIN media, 
and evaluates over piece-wise constant GRIN. One advantage of Trajectron is 
that it is not restricted to the paraxial regime. This makes it eligible for the 
study of non-paraxial optical phenomena such as peripheral visual field analysis, 
as well as higher-order aberrations before and after refractive surgery. 

1.3 Related work 

Previous mathematical models of the human eye assumed the corneal refractive 
index to be constant [21 HH UHl 127] ■ Such models were proposed when less was 
known about the intra-corneal GRIN nature. A few recent models based on 
ray-transfer matrix (ABCD) analysis incorporate the intra-corneal GRIN na- 
ture. Barbero proposed a multilayer model of cornea and tear-film to calculate 
corneal power. His global ray-transfer matrix was a product of 7 matrices rep- 
resenting ray transfer through 4 interfaces and 3 media which included a GRIN 
stroma [5 . Flores- Arias et al proposed a linear cascade model expressing the 
generalized fresnel integral transform kernel in terms of paraxial coefficients; 
where the global ray-transfer matrix was a product of anterior chamber, GRIN 
cornea, and GRIN lens component ray-transfer matrices [m [32] . These newer 
models incorporating the GRIN nature represent a significant advance. One 
disadvantage however, is that they are restricted to the paraxial regime, and 
this limits their validity for peripheral visual field analysis and non-paraxial 
photon tracking. In contrast, our current work is useful for investigating inter- 
esting non-paraxial optical phenomenology, because it respects the intra-corneal 
GRIN nature, and is not restricted to the paraxial regime. 

1.4 Contributions 

In this paper we make the following technical and biophysical contributions: 
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• The first quantitative demonstration of non-paraxial light-bending within 
the human cornea. 

• An analytico-numerical algorithm, Trajectron, for determining the photon 
trajectory within the human cornea, respecting axial GRIN and assuming 
zero meridional GRIN. 

• A method for studying non-paraxial optical phenomena in the human 
cornea. For example: 

— Peripheral visual field analysis 

— Higher order aberrations pre and post refractive surgery 

• An approach towards describing the angular domain over which the parax- 
ial approximation is clinically useful in the human cornea. 

1.5 Methods 

We use a hybrid approach inter-weaving between analytical and numerical pro- 
cesses. Our objective is to determine the ray path from a point in the anterior 
central corneal epithelium to a point in the posterior stroma. The human cornea 
has an axial GRIN of non-constant magnitude |31| , and we approximate it as also 
having zero meridional GRIN. This approximation is made because, although 
meridional intra-corneal GRIN has been demonstrated in- vivo in humans [44j . 
it is significantly less than the axial intra-corneal GRIN [3T]. To achieve our 
objective, we introduce Trajectron, a new algorithm which evaluates a closed 
form solution of the ray equation of constant axial GRIN media, over piece- wise 
constant axial GRIN. The details of the Trajectron algorithm are presented in 
Section ([s]) below. Discrete measurements of corneal refractive indices are ob- 
tained from the literature, and a continuous distribution is generated by curve 
fitting. Trajectron is prototyped in MATLAB® (The MathWorks Inc., Natick, 
MA). 

The remainder of this paper is organized as follows: Section ([2| outlines 
the physical formalism of lagrangian optics and includes a known closed form 
solution of the ray equation of constant axial GRIN media; Section ([s]) de- 
scribes structure, thickness, and refractive index measurements within the hu- 
man cornea; Section Q introduces the Trajectron algorithm; Section ([s]) presents 
the results of our numerical experiments; Section ^ discusses the findings, the 
study limitations, and the clinical angular domain of the paraxial and non- 
paraxial regimes; and Section ([t]) concludes the paper. 

2 Physical Formalism 

Fermat's principle of stationary time asserts that light travels along the path 
that takes the least (or greatest) time. In the context of a classical mechanics 
treatment, Fermat's principle is fully described by the lagrangian of geometric 
optics, given by, 

C = n{x,y,z)^/l + x'^ + y'^, (1) 
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where n is the refractive index distribution, x, y, and z are spatial coordi- 
nates, and Q '■= While x and y take the role of generalized coordinates, z 
describes the optical axis and takes the role of parameter. 

The corresponding Euler-Lagrange equations are given by, 

where Equation ([2| encodes three equations: with q = x, q = y, and q = z. 
Applying Equation ^ above to the optical lagrangian, Equation ([T]), yields 
the ray equation, 

n— — Vn, (3) 



ds \ ds ^ 

which represents the equations of motion for propagation of a light ray 
through a medium of refractive index distribution n{x, y, z). In Equation ([s]), r 
is the spatial position vector [x, y, z] of the ray, V is the gradient operator, and ds 
is the arc length along the ray trajectory, and is given by ds = + + y^. 

In the case of constant refractive index gradient, a, along the optical axis 
direction, and zero refractive index gradient along the other two cardinal direc- 
tions. Equation ^ can be written as. 



d f dr 



ds \^ds 



(4) 



The ray r = [x, ?/, z] is embedded in a 2-dimensional plane. This plane 
of propagation is defined by the optical axis and ds\ ., where i denotes the 
point along the optical axis in the anterior corneal epithelium where we initiate 
observation, and ds\^ denotes the ray direction at i. The embedding, (x, y, z) — > 
{r,z), is given by, 

x — r cos(-0), 

y = rsin(V'), (5) 
z = z, 

where tp is the angle the propagation plane makes with the positive y-axis. 
The components of the canonical optical momentum are given by, 

dx dy dz 

Pi=n--, P2=n--, and = n— , (6) 

as ds ds 

and as follows for direction cosines, sum of squares equals one, hence, 

IpP - n'- (7) 

Equations Q and ^ together state that pi and p2 are constant along the 
ray trajectory. Equation Q in conjunction with the constancy of pi and p2, 
constrain the ray r to the (r, z) plane. 

Applying the fundamental theorem of calculus to the conservation relation 
for optical momentum shown in Equation ([t]), one obtains the following known 
closed form solution of the ray equation in constant axial GRIN media [33j , 
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r(z) = A\n {2a^z + 2aG + 2c 



(8) 



where n.i is the refractive index at the initial observation point of the photon, 
G is a function of z given by, 



G{z) — \J o?z'^ + 2aniZ + nf sin^((/)i), 
A is a constant given by, 

Tli cos(0i) 



A 



(9) 



(10) 



and is the angle between ds]^ and the positive r axis. 

The angle between the optical axis and ds|^ is denoted by Qi. In this paper 
we are primarily interested in rays entering the cornea, and will therefore focus 
on photon trajectories for which < 0^ < | . 



3 Measurement 

Closed-form and numerical solutions of the ray equation of axial GRIN media 
require knowledge of the GRIN. The gradients can be approximated by the 
difference quotients. This requires appropriately sampled measurements of re- 
fractive index and thickness. Below, we examine the literature for the required 
measurements in the four corneal layers along the photon trajectory of inter- 
est: (i) corneal epithelium, (ii) Bowman's layer, (iii) anterior stroma, and (iv) 
posterior stroma. 



3.1 Corneal epithelium 

The corneal epithelium is a non-keratinized stratified squamous epithelium |21j 
of 5-7 cell layers with rapid surface turn over. The corneal epithelium is covered 
by mucinous tear film anteriorly, and is adjacent to Bowman's layer posteriorly. 
The morphology of the air-tear film interface approximately conforms to the 
corneal epithelial surface. Corneal epithelial surface layer abrasions can cause 
significant visual distortion, and can manifest photophobia and moderate to 
severe amounts of pain. However as a result of limbal stem cell activity "26j, 
abrasions typically heal in 24-48 hours in healthy individuals. Both LASEK 
and PRK procedures involve manipulation of the corneal epithelium. In PRK, 
the epithelial layer is removed entirely by mechanical or alcohol debridement, 
to allow for laser ablation of the anterior stroma [MllZl- In LASEK, an alcohol 
solution is used to disrupt the epithelial layer, creating a flap which is then 
replaced following laser ablation of anterior stroma [3l [7j . 



3.1.1 Corneal epithelial thickness 

There have been several measurements of human central corneal epithelial thick- 
ness by various modalities. Significant variation exists in the findings, but al- 
most all the results have been in the 49-60 /im range. The following is a sam- 
ple of prior work in this area. King-Smith et al used interferometry to make 
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measurements on 6 healthy eyes, and found a corneal epithelial thickness of 
49.7 ^m (±2.2 ^m) [2D]; Li et al used Confocal Microscopy Through Focusing 
(CMTF) on 7 healthy subjects and found a thickness of 50.6 //m (±3.9 /im) [25] : 
and Reinstein et al used Very High Frequency Ultrasound (VHFU) on 110 eyes 
in 56 subjects who presented for refractive surgery, and they found a thick- 
ness of 53.4 /xm (±4.6 /xm) [55]. Ladage et al conducted a prospective ran- 
domized double-blind study on 246 subjects, assessing the effects of contact 
lens type and oxygen transmissibility on corneal epithelium in daily wear [23j . 
The study consisted of three treatment groups: high oxygen-transmissible soft 
lenses (n= 36 subjects), hyper oxygen-transmissible soft lenses (n=135 sub- 
jects), and hyper oxygen-transmissible rigid gas permeable lenses (n=75 sub- 
jects), in each of which baseline central corneal epithelial thickness was mea- 
sured using CMTF and found to be 48.55 ijui (±3.2 /zm), 49.26 /xm (±3.44 ^m), 
and 49.79 /im (±3.25 /xm) respectively. There have been a number of studies 
using optical coherence tomography (OCT) [TH |43] to estimate central corneal 
epithelial thickness in humans. Wang et al used OCT on 20 healthy eyes in 
40 subjects and found a thickness of 57.8 /xm (±1.7 /xm) Feng et al used 

OCT on 10 healthy subjects and found a thickness of 58.4 /Ltm (±2.5 /xm); and 
Radhakrishnan et al used 1310 nm real-time OCT to make measurements on 5 
healthy subjects, and found a thickness of 55 /xm [34]. Table summarizes 
the central corneal epithelial thickness for the above mentioned studies. 



Table 1: Thickness of human central corneal epithelium 



Ref 


Method 


Sample size 


Thickness 


M 


Interferomctry 


6 eyes 


49.7 /xm (±2.2 /im) 


131] 


Real-time OCT 


5 subjects 


55 /xm 


[46] 


OCT 


20 eyes 


57.8 /xm (±1.7 /tm) 




OCT 


10 subjects 


58.4 /tm (±2.5 /im) 




CMTF 


7 subjects 


50.6 fim (±3.9 /tm) 


[23] 


CMTF 


36 subjects 


48.55 fim (±3.2 /xm) 


M 


CMTF 


135 subjects 


49.26 /im (±3.44 /im) 


US] 


CMTF 


75 subjects 


49.79 /im (±3.25 /tm) 


m 


VHFU 


56 subjects 


53.4 fim (±4.6 /im) 



3.1.2 Corneal epithelial refractive index 

Patel et al used a modified hand-held refractomcter to measure the refractive 
index in 10 human eyes in-vivo, and found a value of 1.401 (±0.005) [31]. There 
have been few if any other in-vivo estimations of corneal epithelium refractive 
index in humans. 

3.2 Bowman's layer 

In humans. Bowman's layer sits posterior to the corneal epithelium and contains 
randomly arranged collagen fibrils, which on transmission electron microscopy 
(TEM), appear to be continuous with the more regularly arranged fibrils of the 
stroma. Fibrils in Bowman's layer were only half to two thirds the diameter of 
those in the stroma, suggesting that Bowman's layer fibrils may be a condensa- 
tion of stromal fibrils [TH]. Using TEM, Hayashi et al found the collagen fibril 
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diameter to be 20 nm [T7] in Bowman's layer and 20-30 nm in the stroma [T7] : 
and Komai et al found fibril diameters of 25 nm in Bowman's layer and 25-35 
nm in the stroma [22]. 

3.2.1 Bowman's layer thickness 

There have been few measurements of human Bowman's layer thickness, most of 
which have fallen in the 10 — 17 range. In Hayashi et al's TEM studies, mea- 
surements were made on an excised conical cornea and found a central thickness 
of 10 /im. Li et al's CMTF studies on 7 healthy subjects found a thickness of 
16.6 /im (±1.1 /im)d |25j : King-Smith et al's interferometry measurements on 
6 healthy eyes found a thickness of 14.6 /im (±1.4 /im) 20J; and Komai et al's 
TEM studies on 8 eye-bank eyes found a thickness of 8-12 fj,m [22]. Table ^ 
summarizes results from the above sample of work in this area. 



Table 2: Central thickness of Bowman's layer 



Ref 


Method 


Sample size 


Thickness 


M 


Interferometry 


6 eyes 


14.6 fim (±1.4 /im) 


USi 


CMTF 


7 subjects 


16.6 fim (±1.1 /im) 


na 


TEM 


1 eye 


10 /im 


m 


TEM 


8 eyes 


8 - 12 /im 



3.2.2 Bowman's layer refractive index 

There have been no reported measurements of the refractive index function of 
Bowman's layer in humans. However, anterior stroma is intermediate between 
Bowman's layer and posterior stroma, in terms of several structural features such 
as axial depth, susceptibility to hydration swelling [2^ , degree of hydration [23] , 
collagen fibril diameter }17| . and regularity of fibril arrangement [15j . It is 
therefore reasonable to conjecture that the refractive indicex function of the 
anterior stroma is intermediate between those of Bowman's layer and posterior 
stroma. This is premised on the notion that the group refractive index of a 
substance is a function of the refractive indices of the biochemical constituents 
of that substance. Hence the group refractive index can be estimated by a 
weighted summation rule such as the Gladstone-Dale relation; or it can be 
estimated by interpolation, as we do here. 

3.3 Anterior corneal stroma 

The stroma constitutes approximately 90% of the corneal thickness, and the 
anterior stroma is the major target for ablation and remodelling in most refrac- 
tive surgery procedures. The anterior stroma consists of intermediate regularity 
type I collagen fibrils, 20-35 nm in diameter [T7l|22]. On scanning electron mi- 
croscopy (SEM), the collagen fibrils are seen bundled into lamellae which are 
piled in parallel stacks of about 300 lamellae in the central cornea and about 
500 at the limbus [3S]. The anterior stromal lamellae are 0.2-1.2 /im thick and 
0.5-30 /im wide [22], and are interconnected at points either via interlamellae 
crossing collagen fibrils or via anastomosing sublamellae [35], consistent with 
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Maurice and Gallagher's earlier findings in the rabbit cornea [TB]. The archi- 
tecture of the anterior stroma confers shape, curvature, rigidity, and hydration 
resistance to the cornea [29l |6] . 

3.3.1 Anterior corneal stroma thickness 

The anterior corneal stroma has no conventionally defined thickness or demar- 
cation, but the anterior 100 — 120 /im of the stroma has been shown to be 
physiologically distinguishable from more posterior stroma 

3.3.2 Anterior corneal stroma refractive index 

The refractive index of anterior stroma was observed as 1.38 (±0.005) by Patel 
et al, using a bench Abbe refractometer to make measurements on bare stroma 
of fresh human corneas |31j . 

3.4 Posterior corneal stroma 

The posterior stroma differs from anterior stroma in a number of ways. Ul- 
trastructural studies done in rabbit and humans corneas found the posterior 
stroma to have a greater degree of regularity of collagen fibril arrangement than 
the anterior stroma ^15) . Unlike anterior lamellae, posterior lamellae may course 
the entire length of the cornea without branching, and they are typically 1.0 
-2.5 /im thick and 100-200 /im wide, which is larger in dimension than anterior 
lamellae |22j . Horizontal branching occurs among both anterior and posterior 
lamellae, however only anterior lamellae exhibit axial branching [22]. In addi- 
tion, the posterior stroma has been found to be more susceptible to hydration 
and swelling than the anterior stroma [29) . Though the posterior stroma is 
largely spared in current refractive surgery procedures, it becomes the cornea's 
primary source of structural support and bulk following laser ablation of the 
anterior stroma. 

3.5 Posterior stroma thickness 

The posterior stroma can be reasonably designated as the posterior two-thirds to 
three-quarters of the stroma. Reinstein et al used VHFU to measure net stromal 
thickness across the central 10 mm of the cornea in 110 normal eyes, and they 
found a vertex thickness of 465.4 ± 36.9 /im and a thickness of 461.8 ± 37.3 /im 
at the thinnest point [37'. Erie et al used in-vivo confocal microscopy on 18 eyes 
of 12 patients presenting for LASIK, and they found a net stromal thickness of 
491 ± 35 /im pre-operatively [S]. 

3.6 Posterior stroma refractive index 

Patel et al used a bench Abbe refractometer on bare stroma of fresh human 
cornea, and found a posterior stroma refractive index of 1.373 (±0.001). This 
differed from the refractive index of the anterior stroma which they found to be 
1.38 (±0.005) [35. 
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Figure 1: Refractive index sample fit 
3.7 Refractive index distribution 

We fit a curve, shown in Figure Q, to Patel's experimentally measured data 
points [ST shown in Table |3] The fit curve's expression is given by, 

n{z) = 1.3952z-""°3 (n) 

where z > 0.2 is the displacement (in /im) along the optical axis, n{z) is 
the refractive index at z, and we have chosen z = fim as the central corneal 
epithelial surface, z = 0.2 fim as the initial observation point of the photon 
in the corneal epithelium, z = 67 //m as the observation point in the anterior 
stroma, and z — 550 /im as the observation point in the posterior stroma. These 
choices are consistent with the literature data. 



Table 3: Refractive indices in human cornea 



Refractometcr 


Level 


Refractive index 


modified hand-held 
bench Abbe 
bench Abbe 


epithelium 
anterior stroma 
posterior stroma 


1.401 (±0.005) 
1.380 (±0.005) 
1.373 (±0.001) 



4 Computation 

Here we describe the Trajectron algorithm for computing the trajectory of a 
photon in the human cornea, based on geometric optics and the following two 
assumptions: (i) non-zero axial GRIN, and (ii) zero meridional GRIN. 

4.1 The Trajectron Algorithm 

We discretize the optical axis into N segments, with finer segments over Bow- 
man's layer because it has a faster changing refractive index gradient than both 



the corneal epithelium and stroma. To compute the photon trajectory over the 
jth segment we require: (i) the incident refractive index, nj(0), (ii) the incident 
location of the ray, rj(0), (iii) the refractive index difference quotient ctj, and 
(iv) the angle 6j (0) which the incident optical momentum, pj (0) , makes with the 
optical axis. In the above notation we have used the domain parametrization, 
< < 1, over each discrete segment Dj, such that fj denotes the restric- 
tion of an arbitrary function / over the jth segment, fj(0) denotes the value 
of / at the start point of the jth segment, and fj{l) denotes the value of fj at 
endpoint of the jth segment. Domain segment boundary points are related by, 
D^{1) ^ Dj+i{0), for every J. 

(i) We obtain rtj(O) simply by evaluating the fit curve in Equation ( 11 ). (ii) 
The incident ray location is given by. 



r,(0) 



(1), 



(12) 



where rj_i(l), is the ray location at the end of the previous segment and is 
available from the previous computation. 

(iii) The difference quotient is readily computed by, 

_ n,-(l)-n,(0) 

D,il)-D,{Oy ^''^ 

(iv) We then determine the angle the optical momentum p makes with the 
axial direction by. 



where. 



dr 
dz 



j-+i(0) = tan 



1 


dr 






dz 





= A- 



2a- 



+ (a/G) (2ani(0) + 2a2z) 
2a2z + 2ani(0) + 2aG 



(14) 



(15) 



A is a constant defined in Equation ( 10 ), and G is a function of z defined in 
Equation 

It follows that. 



r,{z) = F{j)+r,^,{l), 



(16) 



where F{j) is an evaluation of the closed form solution using the jth seg- 
ment parameters obtained as described above. The Trajectron algorithm is 
summarized in Table Q. 



Table 4: The Trajectron Algorithm 
Input: [D,(0),D,(l)],for j = l,...iV 
Initialize: Ujifi), aj, 6^1 (0), ri(0) = 0, j = 1 

Step 1. F[n,(0),a,,e,(0)] +r,-i(l) ^ 

Step 2. il ^7 

Step 3. tan-^(7) ->• 6ij+i(0) 

Step 4. j ^ i 4- 1 

Step 5. loop 
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5 Results 



We conducted 40 numerical experiments with < 9q,2 < f , each of which 
demonstrated bending of the photon trajectory within the human cornea. Re- 
sults of the first 20 numerical experiments are displayed in Table ([5|. The 
9o.2 values ranged from 2.21 to 79.05 degrees and resulted in ^550 values which 
ranged from 2.25 to 90 degrees. A photon with a 6*0.2 of degrees will proceed 
along a straight line trajectory. However, for all other 60.2 values, some degree 
of bending occurs. Figure ([2| plots both the photon trajectory (solid line) and 
the ray initial tangent (dashed line) for a photon with a 9o,2 of 78.97 degrees. 

12,000 I — ' — • ■ — ' — • ' — ' — ' ' — ' — ' 1 




Corneal z-axis (micrometers) 



Figure 2: Photon trajectory with 6*0.2 of 78.97 degrees and ^550 of 88.67 degrees 
(solid line). Also showing initial tangent direction (dashed line) 

As the 6*0.2 values approach 90 degrees, the bending of the photon trajectory 
increases, and becomes increasingly and exquisitely dependent on angle. The 
difference, i5, between the initial tangent line and the photon trajectory is a 
measure of the extent of bending. We use 6^ to denote the displacement along 
the r-axis at position z /Ltm due to light bending. Table ^ shows the 5z values 
at 70/im, 200/xm, 400/Ltm, and 550/zm. In the next 20 numerical experiments, 
we computed photon trajectories with 60,2 values clustered closer to 90 degrees. 
Specifically, the 80.2 values ranged from 77.68 to 78.72 degrees and resulted in 
6*550 values ranging from 85 degrees to 89.75 degrees. The photon trajectories 
were noted to bend away from the optical axis in a concave upwards direction 
as shown in Figure ([3|. 

6 Discussion 

Until recently, schematic eye models assumed a constant or average corneal re- 
fractive index HTJ [551 HZ] ■ Since emerging evidence of corneal GRIN nature, 
a few newer models have respected the GRIN. But like their older counterparts, 
these newer models are based on paraxial approximation and ray transfer ma- 
trix analysis |14[ I32[ [5]. The first-order paraxial approximation of the wave 
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Corneal z-axis (micrometers) 



Figure 3: Light broom: 20 photon trajectories with 9q,2 ranging from 77.68 to 
78.72 degrees and 6'55o ranging from 85 to 89.75 degrees. 

equation in shortwavelength hmit, is one whose sohitions are rays propagating 
along straight hues. And this hnearity ahows an optical system to be mod- 
elled using ray-transfer matrix analysis. The fitness of these ray-transfer matrix 
analyses depend on the magnitude of high frequency coefficients in the fourier 
transform of the refractive index function. In the case of highly varying re- 
fractive index gradient and consequently highly varying ray spatial gradient, 
linear approximations have limited validity. And a prohibitively large number 
of optical elements or ray transfer matrices may be required to obtain a good 
approximation of the ray path. The transition zone between the corneal ep- 
ithelium and Bowman's layer is an example of such a steep gradient, where 
finer linearization granularity is required. There is not a general notion of when 
linear approximation approaches may begin to lose validity, as the small an- 
gle approximation, sm(9) 0, has no defined boundary. For this reason the 
angular domain over which the paraxial approximation is clinically applicable 
is not clear. Yet clinically, refractive error measurements and correction are 
based on paraxial approximation. Results of our numerical experiments shown 
m Figure ^ may help us gain some traction on this issue. We see an increase 
in non-linearity with 0o.2- This increase in non-linearity is correlated with a loss 
of paraxial approximation validity. Also from Figure ([s]), we see that linear-like 
form is preserved over a much wider angle than only those for which sin(0) « 9. 
This suggests that the angular domain over which the paraxial approximation 
is clinically applicable is more broad than suggested by the small angle approx- 
imation alone. For instance, for a 9q 2 of 50 degrees, the (S550 is 25 /xm, which is 
a relatively small displacement. In addition to exhibiting such preservation of 
linear-like form over broad angle, the paraxial approximation is a potent sim- 
plification device in theory and practice. For these reasons, it has enjoyed much 
success and longevity in the clinic. However, its domain of validity is limited. 
For instance, for 9q 2 > 75 degrees, the corresponding S^^q exceeds 0.5 mm and 
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Table 5: Photon trajectories in the human cornea 



6*0.2 


^550 


570 


^200 


^400 


^550 


(in degrees) 


(in degrees) 


(in /im) 


(in /im) 


(in ^m) 


(in /im) 


79.05 


90.00 


58.24 


299.79 


916.87 


9728.70 


78.84 


87.80 


53.65 


270.30 


801.77 


4983.10 


78.25 


85.71 


43.60 


210.14 


590.23 


2839.90 


76.36 


81.82 


24.82 


110.96 


289.03 


1144.10 


74.00 


78.26 


14.20 


060.90 


153.21 


0567.50 


71.50 


75.00 


08.75 


036.71 


090.79 


0326.45 


66.63 


69.23 


04.12 


016.94 


041.25 


0144.63 


62.20 


64.29 


02.39 


009.71 


023.49 


0081.53 


58.24 


60.00 


01.58 


006.38 


015.38 


0053.10 


54.72 


56.25 


01.14 


004.59 


011.03 


0037.98 


50.14 


51.43 


00.78 


003.13 


007.51 


0025.76 


43.97 


45.00 


00.49 


001.98 


004.74 


0016.25 


35.25 


36.00 


00.28 


001.12 


002.68 


0009.15 


29.40 


30.00 


00.19 


000.78 


001.89 


0006.38 


22.07 


22.50 


00.12 


000.50 


001.18 


0004.04 


17.66 


18.00 


00.09 


000.37 


000.88 


0003.00 


11.78 


12.00 


00.06 


000.23 


000.54 


0001.86 


08.83 


09.00 


00.04 


000.17 


000.40 


0001.36 


04.42 


04.50 


00.02 


000.08 


000.19 


0000.66 


02.21 


02.25 


00.01 


000.04 


000.10 


0000.33 



increases rapidly and non-linearly with 9o.2', such that by 0o.2 = 78.84 degrees, 
the (5550 has reached 5 mm, which is almost the entire radius of the cornea. 
This is clearly non-paraxial territory over which the paraxial approximation is 
invalid. Clinically, this zone is at the fringes of peripheral vision. Loss of vision 
here is not likely to interfere significantly with daily life, and may even go un- 
noticed. However, such subtle peripheral vision loss is often the first sign of an 
insidious optic neuropathy such as primary open angle glaucoma. 

In the numerical experiments presented here, the scope of the photon tra- 
jectories studied are exclusively intra-corneal, with z € [0.2, 550] /im. Based on 
these purely intra-corneal photon trajectories, we have made inference involv- 
ing the visual field. This is justified by a reasonable assumption of continuity. 
We have assumed that there exists a correspondence between the intra-corneal 
angle, 0, and the 9 of the originating ray source in the visual field. In other 
words, relatively wide angles in the field give rise to relatively wide angles within 
the cornea, and relatively narrow angles in the field give rise to relatively nar- 
row angles within the cornea. However, while we have relied on this angular 
correspondence to make inference on field source, it is important to note that 
what we have shown in this paper is that this correspondence is non-linear. 
And more specifically, we have shown that wide intra-corneal angles give rise to 
wider intra-corneal angles, while narrow intra-corneal angles remain relatively 
constant. 

Limitations of the current study include the assumption of gradient zero 
in the meridional directions. Vasudevan et al showed naso-temporal gradients 
in the human corneal epithelium |44) . Using a portable Abbe refractometer 
mounted on a slit-lamp to make measurements on 10 human subjects, they 
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found corneal vertex, nasal, and temporal refractive indices of 1.397 (±0.001), 
1.394 (±0.001), and 1.394 (±0.001) respectively. These gradients are signifi- 
cantly less than those found by Patel [31] in the axial direction, and hence were 
assumed to be zero. The naso-temporal gradients were in the direction of the 
corneal vertex, and their incorporation into the model can be expected to atten- 
uate the concavity of the photon trajectory to some degree. However, the overall 
non-linearity and angular-dependence behavior shown here will be unchanged. 
A second limitation of the current work is the non-uniform specimen prepa- 
ration and measurement method of the axial refractive index. The epithelial 
refractive index was measured in-vivo using a modified hand-held refractome- 
ter, while the stroma was prepared as bare fresh cornea and measured in-vitro 
using a bench-top Abbe refractometer [31] . A third limitation of the current 
study is the sparsity of refractive index data points currently available to us in 
the literature. This sparsity is due to the technological difficulty of making fine 
high precision refractive index measurements in the human cornea. Most of this 
challenging work has been done by Patel and by Vasudevan. Both Patel's ax- 
ial measurements and Vasudevan's naso-temporal measurements contained only 
three data points each, which points towards the challenging nature of the work. 
Unfortunately, mathematical models based on such sparse data may be suscep- 
tible to the aliasing problem. For instance, in the current study we only had one 
data point available for the refractive index gradient in the corneal epithelium. 
However, the corneal epithelium is 55 /zm thick on average with 5-7 cell layers 
and contains histophysiological gradients. The cells are older, flatter, and closer 
to apoptosis anteriorly than posteriorly The epithelium is bathed in the tear 
film mucinous phase anteriorly and continuous with bowman's layer posteriorly. 
These histophysiological gradients suggest that the corneal epithelium contains 
non-zero refractive index gradients. If a gradient does exists, its direction and 
form beg investigation. More generally, finer sampling granularity of the refrac- 
tive index is needed in the human cornea in each of the three principal axis 
directions. In-spite of the above limitations, the current study has provided a 
clear first demonstration of non-paraxial intra-corneal light-bending, and high- 
lighted potential areas for further experimental and theoretical collaboration in 
addressing the global disease burden of refractive error. 

7 Conclusion 

In this paper, we have provided the first quantitative demonstration of non- 
paraxial light-bending within the human cornea. We have also introduced the 
Trajectron algorithm for computing photon trajectories in the human cornea. 
And we have initiated steps towards describing the angular domain over which 
the paraxial approximation is clinically applicable. These demonstrations have 
significant implications in the development of refractive surgery algorithms. 
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